/*

Author : Vishwajeet Dusane
*/

//#include<stdio.h>
//#include<math.h>
#include"Matrix.h"

class Gauss:public Matrix
{


 public:
  double Xn[4]; 
  void step1();      //divide first row by a[1][1] -> Row 1 
  void step2(int ,int );     // subtract R[i][1]-a[i][i]R[1][j]
  void step3(int);     //divide Ri row by a[i][i]
  int  isTransformed();  //returns TRUE if the augmented matrix is Transformed matrix
  void getXnFromTransformedMatrix(); //calculate X1,X2,....,Xn Using
                    		     // Xn=bn
				     // Xn-1=bn-1 - Un-1,n * Xn
				     
  void applyGauss();   // Uses the Guasian Ellimination Algorithm
  
   
};

void Gauss:: getXnFromTransformedMatrix()
{
 
  for(int i=0;i<3;i++)
  {
   Xn[i]=0;
   for(int j=i+1;j<4;j++)
   {
    Xn[i]=Xn[i]+matrixA[i][j];   
   } 
  Xn[i]=matrixB[i]-Xn[i]; 
 }
 
 Xn[3]=matrixB[3];
 
} 
void Gauss::applyGauss()
{

  getMatrixA();
  getMatrixB();
  printMatrix();
  step1();
  
  for(int i=0;i<3;i++)
  {
   for(int j=i+1;j<4;j++)
   {
     step2(j,i); // (forRow,byRow)
   }

   printMatrix(); 
   step3(i+1); 
  } 
      
  if(isTransformed()==1)
   cout<<"Its Working";
  else
   cout<<"Something's Wrong";    
   
 printMatrix();  
 
 getXnFromTransformedMatrix();
 
 for(int i=0;i<4;i++)
  cout<<"X"<<i<<" = "<<Xn[i]<<endl;
}


int Gauss :: isTransformed()
{
    
  for(int i=0;i<4;i++)    //critical loop
  {
    for(int j=0;j<(i-1);j++)
    {
      if(matrixA[i][j]!=0)
      {
      cout<<"In isTransformed :"<<matrixA[i][j];
      return 0;
      }
    }
  
  }
  
  for(int i=0;i<4;i++)    //critical loop
  {
      if(matrixA[i][i]!=1)
      {
      cout<<"In isTransformed For Diagonal :"<<matrixA[i][i];
      return 0;
      }
     
  }
  
  return 1;
}

void Gauss::step1()
{

 double firstElement;
   firstElement=getElement(0,0);
     
   for(int i=0;i<4;i++)
   
      matrixA[0][i]=(double)matrixA[0][i]/firstElement;
    
      matrixB[0]=(double)matrixB[0]/firstElement;
  

}


void Gauss::step2(int forRow,int byRow)
{
   
double position;
double temp;

 position=matrixA[forRow][byRow];
 
  
  for(int j=0;j<4;j++)
   {
        
    temp=(double)matrixA[byRow][j]*position;
    matrixA[forRow][j]=matrixA[forRow][j] - temp;
   
   
  				    
   } 	
   
			    
    matrixB[forRow]=(double)matrixB[forRow] - (position* matrixB[byRow]);
    
    
}



void Gauss::step3(int row)
{
   
   double element;
   //element=getElement(row,row);
     element=matrixA[row][row];
   for(int j=row;j<4;j++)
      matrixA[row][j]=(double)matrixA[row][j]/element;
     
      //cout<<endl<<"step 3"<<endl;
       //printMatrix();
     matrixB[row]=(double)matrixB[row]/element;
  

}



int main()
{
 Gauss G;
 
 G.applyGauss(); 
return 0;
}

 
